Small Numbers of Vortices in Anisotropic Traps 
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We investigate the appearance of vortices and vortex lattices in two-dimensional, anisotropic and 
rotating Bose-Einstein condensates. Once the anisotropy reaches a critical value, the positions of 
the vortex cores in the ground state are no longer given by an Abrikosov lattice geometry, but by a 
linear arrangement. Using a variational approach, we determine the critical stirring frequency for a 
single vortex as well as the equilibrium positions of a small number of vortices. 
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I. INTRODUCTION 

Quantized vortices are one of the hallmarks of superflu- 
idity and the first observation of vortex lattices in liquid 
''He [l| , and more recently in gaseous Bose-Einstein con- 
densates , has lead to a large interest in studying their 
ground states and dynamical properties. Gaseous Bose- 
Einstein condensates are particularly well suited to the 
study of vortices, as their internal and external system 
parameters, such as the intcr-atomic scattering length 
or the density, can be experimentally controlled. This 
allows access to study a large range of vortex states in 
different parameter regimes. 

The structure of a ground state vortex lattice is given 
by the demand to minimise the energy of the system, 
under the condition of fixed angular momentum. For 
large angular momentum this is best done by forming an 
Abrikosov (triangular) lattice of charge one vortices 0,3. 
In different parameter regimes, however, other lattice 
structures have been shown to appear: symmetry pre- 
serving annular arrangements for small condensates [3], 
rectangular lattices in traps with weak, superposed op- 
tical lattices giant vortices in combined harmonic 
and quartic traps [6|, and linear arrangements in atomic 
waveguide structures [1]. 

While vortices are interesting from a fundamental 
point of view, it has recently been pointed out that the 
winding number of a single vortex can be engineered to 
create a topologically stable quantum bit for applications 
in quantum information [^, llo| . Hallwood et al. have 
shown that to create proper macroscopic superposition 
states of angular momentum one needs to identify states 
which are close in energy, strongly coupling, and, at the 
same time, well separated from any other states [ll|. 
While these are difficult conditions to fulfill in gaseous 
Bose-Einstein condensates with short range interactions, 
the promise of long-lived qubits makes the effort worth- 
while. 

A second fundamental building block of quantum in- 
formation processing is the ability to measure and control 
the interaction between two qubits. While interactions 
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between vortices can be studied in Abrikosov lattices, it 
would be of advantage to find a geometry with a lower 
number of nearest neighbours. Two recent works have 
shown that systems where each vortex has only two near- 
est neighbours can be created in atomic waveguide struc- 
ture [8| or in anisotropic, harmonic traps in the limit of 
weak interactions [l^l ■ 

Here we extend the two works above by considering a 
Bose-Einstein condensate deep in the non-linear regime 
(Thomas Fermi limit) and trapped in an anisotropic trap. 
By numerically determining the minimum energy state, 
we find that for relatively small aspect ratios, A ^ 2, and 
moderate rotational frequencies the rotating Abrikosov 
lattice is no longer the ground state of the condensate 
and, instead, a linear crystal of vortices is formed along 
the soft axis of the trap. Taking advantage of the sym- 
metry of the system, we devise an ansatz for a variational 
calculation in order to predict ground state properties of 
the system. We investigate the critical stirring frequen- 
cies needed to move from the ^ = to ^ = 1 state and the 
locations in the lattice for a small number of vortices of 
equal rotational charge. Note that anisotropic traps in 
the strong rotation limit have recently been investigated 
in [13]. 

The paper is organised as follows. In Sec. |TT]we intro- 
duce our variational wavefunction and calculate the en- 
ergy of the system as a function of anisotropy. In Sec. IIIII 
we determine the critical stirring frequencies for local and 
global stability of a single vortex in an anisotropic trap 
and in Sec. IIVI we calculate the ground state structure 
of a vortex crystal with a small number of vortices. In 
Sec. ?? we briefly discuss the dynamics of a linear vortex 
crystal in an anisotropic trap. In Sec. IVII we conclude. 



II. VARIATIONAL CALCULATION 

A. Energy Functional 

We consider a condensate of N atoms with atomic mass 
m and scattering length a^c, trapped in a potential U{y) 
that is rotating at frequency fi. The energy functional 
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for such a system is given by 
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where L is the angular momentum operator, g^D = 
Airh^asc/m is the coupUng constant and / d^r \^\^ = 1. 
Finding the equihbrium states for the above functional in 
ah generahty is a formidable task and beyond our analyt- 
ical as well as numerical capabilities. Wc will, therefore, 
in the following make some assumptions on the symmetry 
of our system. 

First; we will consider rotation around one of the major 
axes only. Therefore, by choosing the angular momen- 
tum operator to be Lz, we can separate the wave func- 
tion into the xy-plane and the z-direction as '^{x, y, z) = 
ip{x,y)(j){z). Since the contributions to the energy func- 
tional from the kinetic, potential and rotational energies 
from are constant, we need only consider the con- 
tributions made to the coupling constant from the third 
dimension. In the following we will therefore make the 
assumption that the wave function in the z-direction is 
in the Thomas-Fermi limit and rescalc the coupling con- 
stant as g2D = gzDi'nnjJz/2'KhY^'^ . This creates an effec- 
tive, two-dimensional Hamiltonian to work with. Speci- 
fying the anisotropy of the harmonic trapping potential 
in the xy-plane by the parameter A, we get the following 
energy functional 



Ng2D 



d^r\ipf 



(2) 



Variational Ansatz 



In order to find a good variational ansatz wc need to 
consider both the effect the angular momentum has on 
the density of the condensate and the contribution to the 
phase of the condensate. A convenient and transparent 
way of doing this is to use the quantum hydrodynamic 
form and split the wavefunction into its modulus and 
phase [14|] 
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Considering the modulus first, we find that, in the 
Thomas- Fermi limit of large particle numbers, a conden- 
sate carrying vortices can be characterised by two length 
scales. The first is the overall condensate size and is given 
by the Thomas-Fermi radius [l^ 
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where fj, is the chemical potential of the system induing 
the vortex lattice. The size of the vortex cores gives the 
second characteristic length scale. For an isotropic trap 
this is given by [l^ 
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and wc will show later that this is still a valid expres- 
sion in the anisotropic case. Since with increasing parti- 
cle number these quantities are inversely proportional to 
each other, we can separate our ansatz with respect to 
these two scales as 




V' = tpTFlpVL 
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The Thomas- Fermi part, ■0tf, describes the background 
cloud and is given by the well known expression 
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The vortices are described by the function V'yLj which 
will only deviate from unity close to the individual vor- 
tex cores and for which we will use a product of tanh 
functions of variable width n = l/j'core [l6j . The full 
ansatz for the wave-function modulus therefore reads 



IV'I = V'TF W tanh [K|r — akR\ 
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and describes n vortices located at d^. Note that the 
vortex positions are scaled in units of the Thomas- Fermi 
radius in the soft direction i?, which lets us restrict our 
values of to between (trap centre) and 1 (condensate 
edge in the soft direction). Wc also define /iq to be the 
chemical potential of a condensate with no vortices 
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To further reduce the number of free variables in the 
energy functional, we will assume that we are in the 
anisotropic limit in which the vortices are forming a lin- 
ear lattice [1, ■ While this assumption at first glance 
looks only reasonable for extremely anisotropic systems 
where the Thomas-Fermi radius is of the order of the 
vortex core diameter, we have found through numerical 
simulations that linear crystals are already present for 
medium anisotropics if the number of vortices is small. 
This can be qualitatively understood by realising that 
the anisotropics lead to extended low density areas in 
the soft direction and an alignment along this axis there- 
fore minimizes the vortex self energies as well as inter- 
actions. In Fig. m^a) we show the standard deviation of 
the vortex core positions from y ~ for a system of 5 
vortices as a function of anisotropy. For an isotropic trap 
the vortices form a triangular arrangement, which with 
increasing values of A becomes squeezed along the tight 
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C. Variational Energies 



Having the defined the ansatz and energy functional in 
the last section, we can now calculate and minimise the 
different energies involved. All our results for single and 
two vortex systems match up with recent work by Castin 
and Dum [l6| in the isotropic limit. 



FIG. 1: Numerical simulations for a *^Rb gas of TV = 10^ 
atoms with a^c = 4.67 x 10~^m and trapping frequency 
LO = 207rs~^. (a) Standard deviation, Oy, of the vortices from 
the X-axis for 1 < A < 2. For A > 1.9 the vortices form a 
linear crystal and Oy = 0. (b,c,d). Ground states for a con- 
densate with parameters as above in a trap with aspect ratio 
of A = 2.5. Each plot spans a range of — < a; < 36/im 

and —Hum < y < Hum. (b) Ground state density for con- 
densate rotating at f2 = 0.7a; showing 6 vortices, (c) Ground 
state density for condensate rotating at Q — 0.65ti; showing 5 
vortices, (d) phase of the cloud shown in (c). The correspon- 
dence between the phase singularities and the density dips is 
indicated. 



direction. For the relative moderate value of A = 1.9 
the standard deviation <7y goes to zero and remains there 
for all higher values of anisotropy. Figs. [T](b)-(d) show 
these numerically determined ground states for an even 
and an odd number of vortices, illustrating this linear ge- 
ometry. Therefore we assume that dk = (afc,0), i.e. the 
vortices sit along the x-axis of the trap. Furthermore, 
due to the symmetry of the system, we can assume that 
for even numbers of vortices, vortex pairs are located at 
(afe,0) and (—0^,0) and for odd numbers of vortices, a 
single vortex will sit at the centre of the trap and the 
remaining vortices will pair off as in the even case. 

Finally, the phase itself can be broken up into the phase 
of the condensate without vortices, 5*0, and the phase 
from each vortex core of charge qk 
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The first part describes the phase of the anisotropic con- 
densate under rotation without vortices and is given by 
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Since in our system the vortices are sitting along a single 
line along which the phase is constant, we can neglect 
this contribution to the phase and take 
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1. Kinetic Energy 



Since we assume our system to be in the Thomas- Fermi 
limit, we will neglect the spatial variation of the back- 
ground condensate cloud when calculating the kinetic 
energy. We can then rewrite tanh^(x) = 1 — sech^(x) 
and take advantage of the fact that sech^(a;) is non- 
zero only within a small distance from the vortex core. 
This allows us to neglect all terms that involve prod- 
ucts of two different vortex cores, i.e. terms of the form 
sech^(xi)sech^(a;j), (i ^ j) and leaves us the following 
expression for the kinetic energy (l6| 
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The first two terms describe contributions from single 
vortices and the last one accounts for the contributions 
from two vortices. As expected from our ansatz, the 
anisotropy is completely contained in the Thomas-Fermi 
parts. The single vortex terms can be straightforwardly 
integrated to give 
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where c = (41n(2) — l)/6 is a constant. The remaining 
term, £'kin2 ; cannot be integrated analytically and we will 
solve it numerically when discussing interacting vortex 
systems below. 
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2.. Potential and Rotational Energies 
The potential and rotational energies are given by 
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and can be integrated to give 
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To be able to calculate the rotational integral we have 
made one further approximation by considering only the 
area traced out by a circle around the vortex whose radius 
is the same size as the trap in the tight direction. This 
means that we only take into account the cylindrically 
symmetric contribution to the phase close to the vortex 
and neglect the asymmetric contributions from the low- 
density edges of the cloud [l3l • 

The total energy of the system can now be expressed 

as 
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where the single vortex self energy is given by 
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and the interaction between two vortices of identical 
charge qk = qk' = ±l(fc 7^ k') can be foimd from 
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trap. Two factors have to be taken into account, namely 
the existence of two different trapping frequencies and the 
fact that the lower one sets an upper limit for rotational 
stability in a harmonic potential. To study the physics 
of a single vortex, we need only consider the vortex self 
energy, W . Wc therefore minimise cq. (j20p first with 
respect to k and find that the size of a vortex core is 
related only to the chemical potential of the cloud at the 
site of the vortex \\m 
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In particular we note that this expression is independent 
of the strength of the anisotropy which justifies our earlier 
assumption that we can separate two length scales in the 
system. Using eq. ([^^ wc can then express the self energy 
as 
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where again the first line corresponds to the kinetic en- 
ergy, the second corresponds to the potential energy, the 
third corresponds to the rotational energy. From this 
expression for the energy, we can derive the two criti- 
cal stirring frequencies, fli and f2g, corresponding to the 
appearance of a local and global energy minimum, re- 
spectively. 

The local minimum corresponds to the point where W 
changes from having a local maximum to a local mini- 
mum at the centre of the trap and is found by setting 
d'^W = at a = 



r.^^^^' 2A3 ^^/Cmo\ 
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where C = For this value of the stirring frequency 

a vortex becomes locally stable in the condensate. 

The global minimum is the frequency for which which 
the energy of the condensate with a vortex at the centre 
has lower energy than the condensate without the vortex, 
i.e. = at a = 
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III. CRITICAL STIRRING FREQUENCY FOR 
A SINGLE VORTEX 



Let us first consider the critical frequency for a single 
vortex state to become the stable ground state of the 



Both quantities are shown in Fig. [5] as a function of 
anisotropy. One can see that both expressions above are 
not simple extensions of the isotropic version, as they do 
not scale linearly with A. This agrees with earlier results 
showing the appearance of irrotational velocity fields in 
rotating anisotropic traps before vortex nucleation [l7| . 
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FIG. 2: Critical frequency for a single vortex as a function of 
aspect ratio, A. Stirring frequencies are defined in terms of 
the trapping frequency in the soft direction, uj. For stirring 
frequencies greater than lu, the cloud becomes unstable. TV — 
1 X 10^, cj = 207rs~^, and a^c = 4.67 x 10~^m is the s-wave 
scattering lenghth for Rubidium-87. 

IV. EQUILIBRIUM LATTICE POSITIONS 

In the following we will calculate the equilibrium po- 
sitions of small numbers of vortices in anisotropic traps 
as a function of anisotropy. While for symmetry rea- 
sons a stable single vortex will sit in the centre of the 
harmonic trap, the positions of larger numbers of vor- 
tices are determined by a complex interplay between the 
single vortex energies of eq. (PO)) . and the vortex- vortex 
interactions given by eq. ([2T|) . As argued above, we will 
assume that the vortices are arranged in a linear crys- 
tal which is aligned along the x-axis and symmetric with 
respect to the y-axis of the trap. 

For two vortices of equal charge and located at ±a the 
integral for the vortex- vortex interaction, eq. (j2ip . can 
be explicitly written out as 

T _ ]^ f ,2 II |2 + - 

~ 2mJ l"^^^' {{x - a)2 + y2) ((j. + a)^ + y^) ' 

(26) 

This expression cannot be fully integrated analytically 
and we therefore minimise the complete energy func- 
tional numerical with respect to the vortex positions ±a. 
The results are shown for different rotational frequencies 
in Fig. [3] as a a function of anisotropy. The regions in 
which the two- vortex molecule is the stable ground state 
of the system are marked by the black parts of the curve. 
As the above minimisation procedure does not give any 
information about the stability of the system, we have 
determined these areas through numerical ground state 
calculations. Note that they are only indicative though, 
since pinning down exact borders is a challenging nu- 
merical task beyond our capabilities. Also note that we 




FIG. 3: Distance of each vortex in a vortex pair from the 
centre of the trap as a function of aspect ratio. The results 
are shown for a gas of TV = 1 x 10® atoms of *^Rb. The 
trapping frequency is given hy uj — 20-ks~^ and we display 
three different rotational frequencies Q = OAu;, 0.6a; and O.Sco. 
The regions in which the two- vortex state is stable is marked 
by the black part of the graphs. 



display real units here, since the Thomas-Fermi radius 
changes as a function of the rotation frequency. 

The first thing one can see from Fig. [3] is that with 
increasing anisotropy the vortices are moving away from 
the trap centre. This behaviour is a reaction to the in- 
creased non-linearities due to the tighter potential in the 
transversal direction, which has the effect of higher sin- 
gle vortex energies as well as vortex-vortex interactions. 
Increasing the rotation frequency has the opposite effect, 
as one can see from the relative position of the three 
curves displayed. Faster rotation, i.e. larger fl, leads to 
an increased Thomas-Fermi radius in the soft direction 
and allows for lower densities near the centre. The vor- 
tices therefore move back towards the centre. The de- 
tailed form of this behaviour for continuous Q is shown 
in Fig. 12 

For the three vortex case we again use the symmetry 
argument to choose an ansatz in which the central vortex 
is located at (0, 0) and the other vortices are located at 
(±a,0). We now have include the interaction between 
all three vortices, giving us three integrals, two of which 
will be identical due to the symmetry. Fig. [5] shows the 
location of the outer vortices as a function of the aspect 
ratio. As with the two vortex case, the vortices settle 
further away from the centre the more anisotropic the 
trap is, however the distance is larger than in the pre- 
vious case. This is clearly due to the addition of the 
repulsive centre vortex and when comparing the two and 
three vortex case, one can see that the distance from the 
centre is about doubled. This indicates that for at least 
small numbers of vortices within the TF cloud the dis- 
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FIG. 4: Distance of each vortex in a vortex pair from the 
centre of the trap as a function of the stirring frequency. The 
system parameters are the same as in Fig. [3] and we display 
four different aspect ratios. The stable regions are again the 
parts marked by the black lines. 

tance between two vortices is almost constant and which 
confirms the overall trend that can be seen from the sim- 
ulations shown in Fig. [T] This behaviour is also known 
from the two-dimensional Abrikosov lattices Q]- Again, 
as in the two vortex case and for the same reason, the 
vortices moves inward with increasing trapping frequency 
(see Fig.E]). 

It is easy to see how this method can be expanded for 
larger numbers of vortices, taking advantage of the sym- 
metry of the system. However, the number of interaction 
terms required for N vortices is N{N — l)/2 requiring 
increased computational power and time. 

V. VORTEX DYNAMICS 

While the above determines the ground state of the 
condensate, it does not guarantee that the dynamical 
evolution maintains the linear shape of the vortex crystal 
over time. In fact, the force on a vortex in an inhomoge- 
neous potential is directed perpendicular to the gradient 
of the potential and one can therefore expect movement 
of each vortex along an cquipotential line of the trap. 
While in an isotropic trap this would lead to a simple ro- 
tation of the linear crystal around the centre of the trap, 
in an anisotropic trap this leads to a change in relative 
distance between two vortices. Shorter distances, how- 
ever, mean larger repulsive forces and one can see that 
beyond a critical anisotropy the vortices will be unable 
to pass each other along the perpendicular direction and 
will therefore maintain the linear crystal structure. 

While the dynamical evolution is not the topic of the 
current work, we have confirmed the above intuition by 




FIG. 5: Distance of outer vortices from central vortex in the 
three vortex case as a function of anisotropy for a range of 
stirring frequencies. The system parameters are the same as 
in the two vortex case. 
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FIG. 6: Distance of outer vortices from central vortex in the 
three vortex case as a function of the stirring frequency for a 
range of anisotropics. The system parameters are the same 
as in the two vortex case. 



simulating the dynamical behaviour of our system for a 
wide range of parameters and anisotropics. In Fig. [7] we 
show the trajectories of five vortex cores in a trap with 
an anisotropy parameter of A = 2.5 for a duration of 
t = 3.5s. One can see that the vortices move around in 
the vicinity of their original position, however they are 
not able to pass each other out. In fact, the linear shape 
of the vortex crystal is maintained at any time and we do 
not even observe any bending. At some point, for even 
longer timescales, the backaction of the vortex movement 
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FIG. 7: Trajectories of the single vortices over a time interval 
of t = 3.5s. One can see that the outer vortices move around 
their equihbrium positions (white dots), but do not pass each 
other out. The central vortex remains stationary throughout. 
The system parameters are the same as in the figures above. 

will lead to excitations of the background density and an 
exchange in vortex position might happen. This will be 
the investigation of a future work. 

VI. CONCLUSION 

Our work fills the gap between two recent investiga- 
tions into vortex lattice structures in anisotropic traps, 
one in the low density limit [l^l and the other in the fast 
rotational limit [l^. We have shown numerically that 
a small change in the anisotropy of a moderately rotat- 
ing anisotropic Bose-Einstein condensate in the Thomas- 
Fermi limit changes the geometry of vortex lattice from 
hexagonal to linear and we have presented a variational 



analysis of such a system. Using the simple symmetry 
of the linear geometry we have devised a straightforward 
ansatz in the Thomas- Fermi limit, which allowed us to 
split the energy function into single vortex terms and two 
vortex interaction terms. From this we have been able to 
determine the critical frequency for the creation of a sin- 
gle vortex and found that as the aspect ratio increases, 
the critical stirring frequency increases. 

Minimising the energy functional for two and three vor- 
tex state, we were able to determine the exact positions 
of the vortex cores. This is a very tedious exercise for di- 
rect numerical treatments, as many vortex lattices config- 
urations have energies closely related to each other. We 
found that the distance of the vortices from the centre in- 
creases with increasing aspect ratio, however, it decreases 
as the stirring frequency increases. At least for systems 
with small numbers of vortices we found that the distance 
between neighbouring cores is effectively constant. 

Finally we briefly addressed the dynamical behaviour 
of a linear vortex crystal and showed that for sufficiently 
anisotropic traps the linear crystal geometry is main- 
tained for very long time-scales. 
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